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Abstract Homicide investigations often depend on the determination of a minimum 
post-mortem interval (PMI m j„) by forensic entomologists. The age of the most devel- 
oped insect larvae (mostly blow fly larvae) gives reasonably reliable information about 
the minimum time a person has been dead. Methods such as isomegalen diagrams or 
ADH calculations can have problems in their reliability, so we established in this study 
a new growth model to calculate the larval age of Lucilia sericata (Meigen 1826). This is 
based on the actual non-linear development of the blow fly and is designed to include 
uncertainties, e.g. for temperature values from the crime scene. We used published 
data for the development of L. sericata to estimate non-linear functions describing the 
temperature dependent behavior of each developmental state. For the new model it is 
most important to determine the progress within one developmental state as correctly 
as possible since this affects the accuracy of the PMI estimation by up to 75%. We 
found that PMI calculations based on one mean temperature value differ by up to 65% 
from PMIs based on an 12-hourly time temperature profile. Differences of 2 °C in the 
estimation of the crime scene temperature result in a deviation in PMI calculation of 
15 - 30%. 

Keywords Forensic entomology ■ Growth modeling • Lucilia sericata ■ Development 
rates • Non-linear model 



1 Introduction 



Insect development rates are applied not only in pest control management but also 



in forensic entomology (Greenberg 19911. Several species of dipteran and coleopteran 



families infest decaying material in order to breed offspring, and this includes the colo- 



nization of animal carcasses as well as dead bodies ( Archer and Elger 2003 Grassberger 
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and Frank 2004 Lane 19751. In homicide investigations, determination of the age of 



larvae feeding on a corpse can indicate a minimal post-mortem interval (PMI I( 



(Goff 19931. This is often important in forensic case work (Benecke 19981 



1.1 Life cycle of blow flies 



The general life cycle of blow flies includes four stages: egg stage, larval stage, pupal 
stage and imago stage ( |Tao| (1927] ). During the larval stage, three instars can be sep- 
arated: 1st, 2nd and 3rd instar, where the latter is divided due to behavioral changes 
in feeding and post-feeding larvae. 

Blow flies deposit egg clutches directly on the food substrate, such as a dead body 
(Smith and Wall 19971, in a position where the eggs are protected and in a moist 



environment. This ensures a food supply for the hatching 1st instar larvae. The first 
three instars each undergo a moult to reach the next developmental stage; the stages 
can be distinguished by the number of respiratory slits at the posterior end of the 
larvae. The third instar stage lasts for longer than the first two and is divided in a 
feeding and a post-feeding phase. The latter is a preparation for pupation. Therefore, 
the larvae leave the food source to find a suitable place for pupation, emptying their 
gut ( Arnott and Turner 2008 1 . About one third of the pre-adult development time is 
spent in the post-feeding larval stage (Greenberg 19911. Then pupation sets in and 



the imago develops within the pupal case till eclosion ( Lowne| |1890[ ). This last stage 



persists for about half of the time of the total development. 

The larval growth rate depends on its body temperature, which is directly influ- 
enced by environmental conditions as ambient temperature and the heat generated 



by maggot aggregations (Slone and Gruner 20071. Also, an important detail for PMI 



determination is that each species has its own temperature dependent growth rate. 



1.2 Methods for PMI determination 



In forensic case work, two different methods are frequently used to calculate a PMI. 
The first uses isomegalen or isomorphen diagrams, by which the lengths or the devel- 
opmental stage of the larvae are combined as a function of time and mean ambient 
temperature in a single diagram (Grassberger and Reiter 2001 1. According to its orig- 
inators, this method is optimal only if the body and therefore the larvae were not 
undergoing fluctuating temperatures, e.g. in an enclosed environment where the tem- 
perature was nearly constant. 

The second method of calculating a PMI estimates the Accumulated Degree Days 
or Hours (ADD or ADH). ADH values represent a certain number of "energy hours" 
that are necessary for the development of insect larvae. The degree day or hour concept 
assumes that the developmental rate is proportional to the temperature within a certain 
species-specific temperature range (overview in (Higley and Haskell 2009)). However, 
the relationship of temperature and development rate (reciprocal of development time) 
is typically curvilinear at high and low temperatures and linear only in between. 

The formula for calculating ADH is given by 



ADH — T ■ (0 — 6>o) 



(1) 
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where T is the development time, O is the ambient temperature, and the minimum 
developmental threshold temperature Oq is a species-specific value, the so called devel- 
opment zero, which is the x-intercept, i.e., an extrapolation of the linear approximation 
of the reciprocal of the developmental time. This value has no biological meaning, it is 



the mathematical consequence of using a linear regression analysis ( Higley and Haskell 
[2009] ). 

One basic condition for using the ADH method is that the ADH value for com- 
pleting a developmental stage stays constant within certain temperature thresholds. 
For example a developmental duration for finishing a certain stage of 14 days at 25° C 
results in 238 ADD when a base temperature of 8°C is assumed. A developmental du- 
ration of 19 days at 21°C results in 231 ADD, both ADD-values are in the same range. 
We analyzed a published data-set for the development of Lucilia sericata (Meigen 1826) 
(Grassbc rger and Reiter||2001[ | and calculated the corresponding ADH values for these 
data. Fig. |l| shows the calculated ADH values for a base temperature of 6*o=8°C (as 
calculated by a linear regression analysis for the used data-set). In the figure we see a 
new effect: for the younger and also shorter developmental phases the ADH values are 
nearly constant over the complete range of temperatures, but for the post-feeding and 
the pupal stages the ADH values are strongly temperature dependent. 

In general, the ADH method seems to give good results only when the larvae of 
interest have been exposed to temperatures similar to those used in generating the 
reference value applied in the PMI calculation (Anderson 20011. Moreover, the tem- 
perature range in which the development rate is actually linear is not wide enough 
to cover all temperatures during a typical summer in Germany (see also examples for 
May/June 2008 in Fig.|3|. Furthermore, neither developmental durations nor base tem- 
peratures for development have been calculated for species originating from Germany. 
The method must therefore be used carefully. 

Furthermore, it is highly problematic that uncertainties for temperature measure- 
ments from a crime scene cannot be taken into account by either of the commonly 
used methods for PMI determination. It is difficult to determine the actual temper- 
ature controlling the larvae at a real crime scene. Since temperature is the variable 
that most influences development, it is crucial to consider it as accurately as possible. 
The standard procedure is to use temperatures of the nearest weather station for the 
desired time frame and correct them by applying a regression starting from tempera- 
tures measured at the crime scene, when taking the larvae as evidence (Archer 20041. 
The corrected values still contain uncertainties that cannot be accounted for by the 
methods currently used for PMI determination. No information exists for either model 
about the quality of the method or the error intervals of the calculated PMIs. 



2 New approach for PMI determination 

We analyzed developmental data for L. sericata at different temperatures and fitted 
an individual exponential function for each developmental stage. Data used as input to 
the model were published by Grassberger and Reiter (2001) and represent the minimal 
time in hours to complete each larval phase (egg stage = stage 0, 1st instar = stage 1, 
2nd instar = stage 2, 3rd instar feeding = stage 3, 3rd instar post-feeding = stage 4 and 
pupal stage = stage 5) until eclosion of the adult blow fly. The used data-set is one of 
the rare sets which covers a lot of temperatures and the resulting growth curve seems 
to represent growth behavior well (see original paper). Unfortunately, Grassberger and 
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Reiter do not give any error values for their measurements, so we assumed an error for 
the developmental times of about 1 hour. These authors used 250 g of raw beef liver 
in plastic jars, and placed 100 eggs on the food substrate. The jars were placed in a 
precision incubator. At each temperature regime the procedure was repeated 10 times. 
Every 4 hours, four of the most developed maggots were removed from the plastic jars, 
killed in boiling water, and preserved in alcohol ( jAdams and H all , 2003 ) and then their 
stage of development was determined. 



2.1 Data fit 



Our new larval growth model is based on the data shown in Fig. [2] in which the duration 



of each developmental stage was measured as a function of temperature ( Grassberger 



|and Reiter| |2001[ ). These data points were fitted with an exponential function of the 
form 



T a (9) = a a ■ exp (-t q ■ 9) + T , Q (2) 

where T a is the duration of one developmental stage a as a function of temperature 
9. The parameters fitted for the different stages are shown in Table [TJ The param- 
eter r Q defines how strongly the time interval depends on temperature; the higher 
the parameter in Table [TJ the steeper is the gradient of the fitted curve. Tq^ repre- 
sents the minimum time interval required for finishing a certain developmental stage 
and a a provides the absolute normalization. The developmental stages of the mag- 
gots were determined every AT = 4h, such that time measurement errors are set to 
ot> — Z\T/vT2 following an uniform distribution. It is assumed that the maggot body 
temperature is known to an accuracy of 3 % in order to take into account uncertainties 
about differences between ambient and maggot body temperature. The parameters 
a a , r a and To jQ were determined by minimizing the sum of error squares. As seen 
in Fig. [2] the exponential function accurately models the behavior during all devel- 
opmental stages and will be used below. In all stages the developmental duration at 
temperatures below 24 °C starts to rise exponentially. Fig.[l]shows the calculated ADH 
values corresponding to eq. |IJ) (data points). In addition, the figure shows the function 
ADH a (<9) = T a {9) ■ (6>-6> ) (lines). T a (9) is calculated by eq. (2} with the previously 
fitted parameters (Table [TJ. Again, the functions give a reasonable description of the 
data. Nevertheless, the model is an empirical one, based on the observations of the 
data points generated by Grassberger and Reiter (2001). 



2.2 PMI calculation 



For European and especially German temperatures, calculation of the total develop- 
mental duration must allow for non-linear temperature behavior in order to ensure 
accuracy. The basic idea underlying a new approach in PMI determination is to fol- 
low an ambient time-temperature profile 9(t) backwards in time starting from the 
time point tp at which the maggots of interest were collected. The idea of backwards 
calculation is obviously similar to the ADH method, but in the new model the im- 
portant improvement is the way of calculating the larval age. The latter is calculated 
successively during certain time steps using the fitted functions (introduced in Fig. [2J) 
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corresponding to the current developmental stage. In each stage a the relative devel- 
opmental progress is P a (values 0-1) where is the beginning and 1 is the finishing 
point of each developmental stage; e.g. a maggot in the middle of the post-feeding 
stage is P4 — 0.5, at the end of the post-feeding stage it is P4 = 0.9 and so forth. The 
developmental duration t a fi spent in each individual stage is calculated by solving the 
relation: 

where dt/T(0(t)) is the infinitesimal relative development. 

The calculation starts with the developmental stage of the maggot at the time 
of collection, summing the developmental progress of each stage backwards until the 
beginning of the egg stage is reached. The calculation for each collection stage uses 
£q+i,0 = %F- The total development time to or post-mortem interval (PMI) is then 
given by 

io = ^*a,o- (4) 

a 

For the new model a program was written in C++ using Root (http://root.cern.ch/ll. 
This program includes all mentioned mathematical steps and produces the figures 
shown here as output. For each new PMI calculation, the corresponding temperature 
profile can be inserted and individually chosen uncertainties can be included. 



2.3 Consideration of uncertainties by Monte-Carlo simulation 

To explore the uncertainties in the total developmental duration, a Monte-Carlo sim- 



ulation was applied, which is commonly used for simulations in life sciences (Mansson 



et al. 20051. It is a method for calculating one final uncertainty after considering all 
statistically independent uncertainties that influence e.g. the larval age. The mean 
PMI with corresponding standard variation is calculated n times taking into account 
and varying all uncertainties described in the following. First, the developmental pro- 
files T a (&) have uncertainties due to the measurement procedure. Second, the time- 
temperature profile from the collection scene is not known precisely and must be ap- 
proximated using temperature values from nearby weather stations. The variations are 
introduced for each model as follows: 

Development profile: The mean duration values of the temperature-time data are ran- 
domly smeared with a uniform distribution with corresponding error a?; for the 
maggot body temperature oq ^ a Gaussian distribution is used. New fits with the 
function in eq. (J2J are performed for each stage. 

Time-temperature profile: Deviations between the temperature profile at the collection 
scene and the nearest weather station are accounted for by Gaussian smearing of 
time t and temperature 0, with the corresponding errors at and oq as width 
for each data point. o@ can be inserted in the model's calculation individually 
dependent on the differences between the temperatures at finding place and weather 
station. 

We calculated the PMI for a mock crime scene with the following parameters: the 
error of the measurement of the original data aj< — 4/\/l2h, the errors of the data of 
the weather station at = 1 h and a@ = 2°C, the difference of the ambient temperature 
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and larval body temperature ctqj, — 3% for 10.000 models for a fixed collection stage 
progress of P a = 0.5. The results are shown in Fig. [3] (upper part) which is a direct 
output of the new program that calculates the PMI. The lower time axis defines the 
progress of the temperature profile forward in time, representing the time frame of 
interest. The temperature profile used here (black line) is taken from the minimum and 
maximum temperatures in May and June 2008 measured at Cologne/Bonn airport. The 
right end of the diagram marks a fictional time point of maggot collection and therefore 
the starting point for PMI calculation. The upper time axis depicts the PMI backwards 
in time starting from the moment of maggot collection. For each developmental stage 
the PMI was calculated by following a linear interpolation between the maximum and 
minimum temperatures. The histograms illustrate the PMI distribution for each stage 
and show a clear single peak structure. The arrows on the top show the 1-standard 
deviation interval for each stage around the mean PMI value, and range between 0.1 
and 1.2 days (depending on the stage). Since no data points below temperatures < 
15°C were measured, the functions T a (0) were extrapolated to lower temperatures. 
As expected, the PMI and the corresponding standard deviation increase with higher 
developmental duration (see arrows above histogram). 

Since the exact progress within the developmental stage at collection time is most 
of the time also unknown, a third uncertainty is introduced: 

Stage progress: The developmental stage at collection time was determined only to 
integer precision, so that it is assumed the exact progress is an uniformly distributed 
value between and 1. Consequently, the starting value for the PMI calculation 
P a at time tp is randomly and uniformly chosen within the interval [0,1] for each 
model. 

Fig. [3] (lower part) shows the PMI calculation for the same parameters as before, 
but without setting the progress of the development for each stage to a fixed value. 
The 1-standard deviation values increase by 0.3 to 3.3 days. The resulting uncertainty 
in the progress of the stage contributes about 75 % to the total PMI error interval. 
In addition, the histograms show deviations from a clear single peak structure, e.g. 
for the pupal stage, implying that the PMI probabilities for 21 days and 26 days are 
nearly the same. To use the new model, the crucial parameter is therefore the correct 
determination of the progress of the developmental stage of maggots collected from a 
corpse. 

2.4 Estimation of temperature at the location of maggot collection 

The impact of correct temperature determination at the maggot collection scene is 
shown in Fig. [4] The data points represent the mean PMIs with an error bar of 1 
standard deviation as a function of collection stage for three different temperature 
profiles. The triangles show the PMIs for the original temperature profile as measured 
at Cologne/Bonn airport. The bullets (squares) show the results for the same profile but 
subtracted (and added) by 2 °C. As expected, the PMIs and the corresponding standard 
deviations of the lower (higher) temperature profile increase (decrease) relative to the 
nominal profile. These differences in temperature of 2 °C give rise to an effect of 15 - 
30 %. That implies that a miscalculation of the temperature at the crime scene of 2 °C 
will result in a miscalculation of the PMI by 15-30%. The later the stage, the greater 
the deviation from the actual PMI. 
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2.5 Comparing PMIs based on mean temperatures and 12-hourly temperature profiles 

In Fig. [5j PMIs calculated using the corresponding mean temperature values in the 
temperature interval [tojip] are compared with PMIs from our model using the three 
temperature profiles introduced previously. The calculated mean temperatures were 
as follows (calculated for the time frame till completion of each stage): stage = 
f8°C, stage I = 19 °C, stage 2 = 20 °C, stage 3 = 19 °C, stage 4 = 16 °C, stage 
5 = 17 °C. The PMI values based on the temperature profile and those based on a 
mean temperature value agree to within about 5 % for the high temperature value 
(original profile + 2 °C) in all stages. The deviation between mean temperature and 
the original temperature profile exceeds the 10 % level starting at the 3rd instar feeding 
stage, and increases to 25 % in the pupal stage. This effect becomes even larger for the 
low temperature profile (original profile subtracted by 2 °C). Starting from the 2nd 
instar stage, the deviation increases from about 10 % up to about 65 % for the pupal 
stage. This means that use of mean temperature values overestimates the influence of 
low temperatures and underestimates periods of high temperatures. The effect should 
be larger if the mean temperature during the development is lower still, e.g. in spring 
or fall. In general, more data points are needed for the developmental duration at low 
temperature ranges to provide more reliable statements. 



2.6 Does the model work in a real case? 

We calculated the PMI in a real case where the actual PMI was known due to a 
confession of the offender. At the end of August 2007 the victim was killed in early 
morning and was found 4 days later also in the morning on a grassland. This leads to 
a PMI of approximately 96 hours. The victim was stabbed to death and had several 
wounds which would act as attractant to the blow flies. It can be assumed that blow 



flies started ovipositing early after death occurred (Reibe and Madea 20101. Autopsy 
was performed directly after the corpse was recovered and several 2nd instar larvae of 
L. sericata were collected. The largest larvae measured 6.1 mm. Hourly temperature 
values were taken from a weather station 10 km away. The mean temperature was 16 
°C. Using Grassberger and Reiter's isomegalen diagram for a larvae measuring 6 mm 
and a mean temperature of 16 °C results in a time interval of 3.2 days plus 30 hours 
(larval development time plus egg period). In total a PMI of 107 hours is indicated. 
This would shift the time of oviposition to nighttime, which is a highly unlikely event 
(Amendt et al.l [2008b. The same data can be used to calculate the ADH value for 



L. sericata for reaching 6 mm in order to calculate the PMI not based on the mean 
temperature but on hourly data. As mentioned earlier, a regression analysis of the data 
set reveals a base temperature of 8 °C. The corresponding ADH value is therefore 856, 
based on the equation: 

ADH = 107(16- 8) (5) 

Subtracting the hourly ADH values, estimated by the temperature values from the 
weather station and the base temperature, from the starting value of 856 results in a 
PMI of 101 hours. 

To use the new model for calculating larval age in the real case, information about 
the progress of the 2nd instar larval stage was required. In the original work of Grass- 
berger and Reiter (2001, Figure 1) a figure is included showing the growth of the larvae 
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and also the time points for each moult. According to this figure, the 2nd instar stage 
sets in after the larvae have reached a size of approximately 4 mm and ends when the 
larvae have reached a size of approximately 8mm. As the largest larvae we collected 
measured 6 mm, we chose P=0.5 as progress for the larval stage. We included the 
hourly temperature profile and chose a temperature error of 1 °C. The result of the 
calculation was a PMI of 99 hours (SD = 3 hours). 

These calculations of a PMI in a real case show that all three methods give rea- 
sonable results. Furthermore, it becomes obvious that the new model is a possible 
alternative for the existing methods with the benefit of directly providing a standard 
deviation for the calculation. 



3 Conclusion and outlook 

The new model improves the larval age calculation in specific ways. It can be used 
in non-linear parts of the temperature dependent development, and includes individu- 
ally defined uncertainties for a temperature profile determined retrospectively from the 
nearest weather station. In the new model the temperature profile plus the determina- 
tion of the larval stage are translated into a mean PMI as well as a standard deviation. 
PMI calculation using mean temperatures, however, can lead to severe deviations from 
the real PMI. 

So far, the main uncertainty arises from the fact that the developmental stage is 
determined only on a 1 - 6 scale (egg, 1st instar, 2nd instar, 3rd instar feeding, 3rd instar 
post-feeding and pupae) . As shown above, 75 % of the uncertainties in the model depend 
on the exact determination of the developmental progress, and additional length values, 
as shown for the PMI calculation in the real case, will propably increase its accuracy 
leading to more accurate PMI calculations. Moreover, the next step is to produce own 
growth data with known error values to refine the inclusion of uncertainties that are 
only rough estimates at the present time and to improve the till now only empirical 
model. 

Nevertheless, the new PMI calculation program is suitable for use in forensic case 
work as a general tool for PMI determination. Scientists from every country or climatic 
region can incorporate their own growth values for different species and ensure a high 
accuracy in PMI determination. 
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Fig. 1 Calculated ADH values for the development of L. sericata using eq. Ilk (data points) 
and fitted functions (lines) calculated using eq. (|2j and estimated parametres (Table [TJ . Orig- 
inal data by Grassberger and Reiter 2001. 



Table 1 Fitted parameters for the development-time- function T a (0). 
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Fig. 2 Developmental data of L. sericata with fitted functions (eq. |2|). Original data by 
Grassberger and Reiter 2001. 



12 



♦MM 

PMI [days] 

-40 -35 -30 -25 -20 -15 -10 -5 




time [h] 



■4 — ► 

PMI [days] 

-40 -35 -30 -25 -20 -15 -10 -5 




time [h] 

Fig. 3 Excerpt of the temperature profile during May and June 2008 in Cologne/Bonn, Ger- 
many (www.wetteronline.de). The histograms illustrate the PMI distributions of the different 
random models for a certain collection stage a. The arrows at top show 1-standard deviation 
intervals around the mean PMIs for each stage. Upper part: PMI calculation with parameters 
T = 4/%/l2h, at = 1 h, crg> = 2°C, ff©,6 = 3% for 10.000 models for a fixed collection stage 
progress of P a = 0.5. Lower part: PMI calculation as before but the progression of the stage 
is not fixed. 
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Fig. 4 Calculated PMIs and SDs for each stage for a temperature profile of May/ June 2008 
measured at Cologne/Bonn airport in comparison to profiles zfc 2 °C. The later the stage the 
greater the SD and the differences in the calculated PM1 for differences of 2 °C. 
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Fig. 5 Deviations of calculated PMIs when using the mean temperature and the actual 12- 
hourly temperature profile of May/June 2008 from Cologne/Bonn airport. Additionally, the 
results are shown for the temperature profiles ± 2 °C. 



